% testMonitorTiming
%
% play with timing using PMD-1208FS, Optoelectronics photodiode, a bottle
% of wine and a smile
%

%% flash for one second
SCREEN_INDEX  = 1;  % 0=small rectangle on main screen; 1=main screen; 2=secondary
NUM_REPS      = 10;
NUM_LUMINANCE = 5;
GRAPHICS_TYPE = 'snow_dots'; % 'snow_dots', 'mgl', or 'dotsX'
SHOW_RAW      = true;
SHOW_FITS     = true;
REPORT_FITS   = true;

%% setup screen
if strcmp(GRAPHICS_TYPE, 'snow_dots')
    dotsTheScreen.reset();
    s = dotsTheScreen.theObject;
    s.displayIndex = SCREEN_INDEX;
    dotsTheScreen.openWindow();
    
    % make target on the center of the screen ... flash a few times
    t = dotsDrawableTargets();
    t.xCenter = 0;
    t.yCenter = 0;
    t.width   = 10;
    t.height  = 10;
    for ii = 1:10
        t.colors = ones(1,3).*mod(ii,2);
        dotsDrawable.drawFrame({t});
    end
    
elseif strcmp(GRAPHICS_TYPE, 'mgl')
    
    mglOpen(2);
    mglVisualAngleCoordinates(57,[40 30]);
    frameRate  = ceil(mglGetParam('frameRate'));
    deviceRect = mglGetParam('deviceRect');
    mglClearScreen([0 0 0]);
    mglFlush;

else % dotsX
    
end

%% configure device
aIn = AInScan1208FS();
aIn.channels  = 0; % differential channel 0 is inputs 0 & 1
aIn.gains     = 7; % 20x = +/-1V; see MCCFormatReport
aIn.frequency = 2000;

%% Check dark->light and light-dark transitions
if NUM_LUMINANCE == 0
    lums = 0;
else
    lums = linspace(0,0.4,NUM_LUMINANCE);
end
preFrames     = 8;
postFrames    = 8;
nFrames       = preFrames+postFrames;
nSamples      = ceil((preFrames./ceil(s.windowFrameRate)+0.1)*aIn.frequency);
aIn.nSamples  = nSamples.*4;
frameTiming   = nans(nFrames, NUM_LUMINANCE, 2, NUM_REPS);
OSOdata       = nans(nSamples, 2, NUM_LUMINANCE, 2, NUM_REPS);
for ll = 1:NUM_LUMINANCE
    lum_ld = [lums(ll) 1-lums(ll); 1-lums(ll) lums(ll)];
    for dd = 1:2
        for rr = 1:NUM_REPS
            frames = cat(1,repmat(lum_ld(dd,1),preFrames,1),repmat(lum_ld(dd,2),postFrames,1));
            %s.resetFlushGauge();
            configTime = aIn.prepareToScan();
            startTime  = aIn.startScan(); % host CPU time, when start ack'd by USB device
            for ii = 1:nFrames
                t.colors = ones(1,3).*frames(ii);
                timing   = dotsDrawable.drawFrame({t});
                if ii > 1 && ~timing.isTight
                    disp([ll dd rr ii])
                end
                frameTiming(ii,ll,dd,rr) = timing.onsetTime;
                % mexHID('check');
            end
            pause(0.2);
            mexHID('check'); % probably not necessary...
            
            % cleanup
            stopTime = aIn.stopScan();
            [chans, volts, times, uints] = aIn.getScanWaveform();
            t0 = frameTiming(preFrames+1,ll,dd,rr);
            ti = find(times>=(t0-0.01), 1);
            OSOdata(:, :, ll, dd, rr) = [ ...
                (times(ti+(1:nSamples))'-t0).*1000, ...
                volts(ti+(1:nSamples))'];
            frameTiming(:,ll,dd,rr) = frameTiming(:,ll,dd,rr) - t0; 
        end
    end
end

% clean up
dotsTheScreen.closeWindow();
aIn.close();

%% Plot raw data
if SHOW_RAW
    figure
    hold on
    for ll = 1:NUM_LUMINANCE
        for dd = 1:2
            subplot(NUM_LUMINANCE,2,(ll-1)*2+dd); cla reset; hold on;
            title(sprintf('Luminance = %.2f', lums(ll)))
            for rr = 1:NUM_REPS
                startTime = frameTiming(5,ll,dd,rr);
                plot(OSOdata(:,1,ll,dd,rr), OSOdata(:,2,ll,dd,rr), '-', ...
                    'Color', lums(ll).*ones(1,3))
                for ff = 1:size(frameTiming,1)
                    plot((frameTiming(ff,ll,dd,rr)).*1000.*[1 1], [0 1], 'r-');
                end
                plot((frameTiming(preFrames+1,ll,dd,rr)).*1000.*[1 1], [0 1], 'b-');
            end
            axis([-20 150 0.20 0.36]);
        end
    end
end

%% fit to piecewise function, with exponential
if SHOW_FITS || REPORT_FITS
    
    % get fits
    global data
    inits = [40 0.3 10 0.35; 0 0 0.1 0; 100 1 20 1];
    fits = nans(NUM_REPS,size(inits,2),NUM_LUMINANCE,2);
    for ll = 1:NUM_LUMINANCE
        for rr = 1:NUM_REPS
            for dd=1:2
                Lt = OSOdata(:,1,ll,dd,rr)>0 & OSOdata(:,1,ll,dd,rr)<120;
                times = OSOdata(Lt,1,ll,dd,rr);
                volts = OSOdata(Lt,2,ll,dd,rr);
                yStartg = volts(1);
                yEndg   = volts(end);
                [Y,I]   = max(abs(diff(nanrunmean(volts,5))));
                tStartg = min(times(I));
                data    = OSOdata(Lt,:,ll,dd,rr);
                [ff, fval, exitflag, output] = patternsearch(@expFitErr,...
                    [tStartg yStartg 3 yEndg],[],[],[],[],inits(2,:),inits(3,:),[], ...
                    psoptimset('Display', 'off'));
                fits(rr,:,ll,dd) = ff;
                
                %             cla reset; hold on;
                %             plot(times, volts, 'k-');
                %             plot(times, valFIT_piecewiseEXP(ff(1), ff(2), ff(3), ff(4), times), 'r-');
                %             axis([0 120 0.3 0.38]);
                % r = input('next')
            end
        end
    end

    % report
    if REPORT_FITS
        for ll = 1:NUM_LUMINANCE
            disp(sprintf('Lum=%.2f: dark->light = %.2f ms delay, %.2f rise time', ...
                lums(ll), median(fits(:,1,ll,1)), median(fits(:,3,ll,1))))
            disp(sprintf('Lum=%.2f: light->dark = %.2f ms delay, %.2f rise time', ...
                lums(ll), median(fits(:,1,ll,2)), median(fits(:,3,ll,2))))
        end
    end
    
    if SHOW_FITS
        %% plot fits
        figure
        hold on
        lm = [0 50; 0 20];
        params = [1 3];
        nParams = length(params);
        for pp = 1:nParams
            for dd = 1:2 % dark->light, light->dark
                subplot(nParams,2,(dd-1)*2+pp); cla reset; hold on;
                plot(fits(:,params(pp),1,dd), 'k.', 'MarkerSize', 10);
                %         for ll = 1:NUM_LUMINANCE
                %             plot(ll, fits(:,params(pp),ll,dd), 'k.');
                %             plot(ll, median(fits(:,params(pp),ll,dd)), 'r.', 'MarkerSize', 10);
                %         end
                axis([0 NUM_LUMINANCE+1 lm(pp,:)])
            end
        end
    end
end









%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% DONE


%% Plot aligned curves
clf; hold on;
ti = 1;
vs = reshape(OSOdata(:,2,:,:),[],1);
minv = min(vs)-0.01;
maxv = max(vs)+0.01;
for ll = 1:NUM_LUMINANCE
    for dd = 1:2
        subplot(NUM_LUMINANCE,2,(ll-1)*2+dd); cla reset; hold on;
        plot([0 0], [0 1], 'k:');
        for rr = 1:NUM_REPS
            startTime = min(fits(:,1,ll,dd));
            plot(fits(rr,1,ll,dd)-OSOdata(:,1,ll,dd,rr)+startTime, OSOdata(:,2,ll,dd,rr), '-', ...
                'Color', lums(ll).*ones(1,3))
        end
        axis([-20 100 minv maxv]);
    end
end




%% Frame timing data
% frameTiming   = nans(nFrames,  3, NUM_LUMINANCE, 2, NUM_REPS);
for ll = 1:NUM_LUMINANCE
    for dd = 1:2
        subplot(NUM_LUMINANCE,2,(ll-1)*2+dd); cla reset; hold on;
        plot(squeeze(frameTiming(9,1,ll,dd,:)-frameTiming(1,1,ll,dd,:)), 'r.');
        plot(squeeze(frameTiming(10,2,ll,dd,:)-frameTiming(1,1,ll,dd,:)), 'g.');
        plot(squeeze(frameTiming(10,3,ll,dd,:)-frameTiming(1,1,ll,dd,:)), 'b.');
        ylim([0.1 0.2]);
    end
end


for ll = 1:NUM_LUMINANCE
    for dd = 1:2
        subplot(NUM_LUMINANCE,2,(ll-1)*2+dd); cla reset; hold on;
        for rr = 1:NUM_REPS
            plot(frameTiming(2:end,3,ll,dd,rr)-frameTiming(1:end-1,1,ll,dd,rr), 'r.');
%             plot(frameTiming(:,1,ll,dd,rr)-frameTiming(1,1,ll,dd,rr), 'r.');
%             plot(frameTiming(:,2,ll,dd,rr)-frameTiming(1,1,ll,dd,rr), 'g.');
%             plot(frameTiming(:,3,ll,dd,rr)-frameTiming(1,1,ll,dd,rr), 'b.');
        end
%        axis([0 16 0 0.02])
    end
end

    


%%%%%% OLD STUFF BELOW
%%%%%%

% get stats
% onset delay
delays = nans(NUM_LUMINANCE, NUM_REPS, 2, 2);
for ll = 1:NUM_LUMINANCE
    for rr = 1:NUM_REPS
        for dd=1:2
            ts   = OSOdata(:,1,ll,dd,rr)>0 & OSOdata(:,2,ll,dd,rr)>100;
            minv = min(OSOdata(ts,2,ll,dd,rr));
            maxv = max(OSOdata(ts,2,ll,dd,rr));
            Fonset = find(ts & OSOdata(:,2,ll,dd,rr) >= minv+0.75*(maxv-minv), 1);
            delays(ll,rr,dd,:) = [OSOdata(Fonset,1,ll,dd,rr) OSOdata(Fonset,2,ll,rr)];
        end
    end
end


%% loop through each luminance interval
NUM_LUMINANCE = 10;
lums          = linspace(0,0.45,NUM_LUMINANCE);
fseq          = [1 1 0 0 1 1 1 0 0 0 1 0 1 0];
nFrames       = length(fseq);
nSamples      = ceil((nFrames./ceil(s.windowFrameRate)+0.1)*aIn.frequency);
aIn.nSamples  = nSamples.*3;
frameTiming   = nans(nFrames, 3, NUM_LUMINANCE, NUM_REPS);
OSOdata       = nans(nSamples, 2, NUM_LUMINANCE, NUM_REPS);
for rr = 1:NUM_REPS
    for ll = 1:NUM_LUMINANCE
        % make luminance frame  sequence
        lseq = lums(ll).*ones(size(fseq));
        lseq(~fseq) = 1-lums(ll);
        configTime  = aIn.prepareToScan();
        startTime   = aIn.startScan(); % host CPU time, when start ack'd by USB device
        t.colors = zeros(1,3);
        for ii = 1:10
            timing = dotsDrawable.drawFrame({t});
        end
        for ii = 1:length(lseq)
            t.colors = ones(1,3).*lseq(ii);
            timing   = dotsDrawable.drawFrame({t});
            frameTiming(ii,:,ll,rr) = [...
                timing.onsetTime, ...
                timing.swapTime, ...
                mglGetSecs()];
            mexHID('check');
        end
        pause(0.2);
        mexHID('check');
        
        % cleanup
        stopTime = aIn.stopScan();
        [chans, volts, times, uints] = aIn.getScanWaveform();
        t0 = frameTiming(1,1,ll,rr);
        ti = find(times>=(t0-0.01), 1);
        OSOdata(:, :, ll, rr) = [ ...
            volts(ti+(1:nSamples))', ...
            times(ti+(1:nSamples))'-t0];
    end
end
dotsTheScreen.closeWindow();
closeStatus = aIn.close();

% get stats
% onset delay
delays = nans(NUM_LUMINANCE, NUM_REPS, 2);
for ll = 1:NUM_LUMINANCE
    for rr = 1:NUM_REPS
        minv = min(OSOdata(:,1,ll,rr));
        maxv = max(OSOdata(:,1,ll,rr));
        Fonset         = find(OSOdata(:,1,ll,rr) >= minv+0.75*(maxv-minv), 1);
        delays(ll,rr,:) = [OSOdata(Fonset,2,ll,rr).*1000 OSOdata(Fonset,1,ll,rr)];
    end
end

% plotz
cla reset; hold on;
for ll = 1:NUM_LUMINANCE
    subplot(NUM_LUMINANCE,1,ll); cla reset; hold on;
    for rr = 1:NUM_REPS
        startTime = frameTiming(1,1,ll,rr);
        plot(OSOdata(:,2,ll,rr).*1000, OSOdata(:,1,ll,rr), '-', ...
            'Color', lums(ll).*ones(1,3))
        plot(delays(ll,rr,1), delays(ll,rr,2), 'gx');
        for ff = 1:size(frameTiming,1)
            plot((frameTiming(ff,1,ll,rr)-startTime).*1000.*[1 1], [0 1], 'r-');
        end
    end
    axis([-20 300 0.28 0.36]);
end




cla reset; hold on;
for ll = 1:NUM_LUMINANCE
    for ff = 1:size(frameTiming,1)
        startTime = frameTiming(ff,1,ll);
        plot((OSOdata(:,2,ll)-startTime).*1000, OSOdata(:,1,ll), '-', ...
            'Color', lums(ll).*ones(1,3))
    end
end
%axis([0 100 0.3 0.36])



cla reset; hold on;
for ll = 1:NUM_LUMINANCE
    startTime = frameTiming(1,1,ll);
    plot((OSOdata(:,2,ll)-startTime).*1000, OSOdata(:,1,ll), '-', ...
        'Color', lums(ll).*ones(1,3))
    for ff = 1:size(frameTiming,1)
        plot((frameTiming(ff,1,ll)-startTime).*1000.*[1 1], [0 1], 'r-');
    end
end
ylim([0.3 0.36]);



plot((times-startTime).*1000, volts, 'k-')
for ff = 1:size(frameTiming, 1)
    plot((frameTiming(ff,1)-startTime).*1000.*[1 1], [0 1], 'r-');
    plot((frameTiming(ff,2)-startTime).*1000.*[1 1], [0 1], 'g--');
    plot((frameTiming(ff,3)-startTime).*1000.*[1 1], [0 1], 'b--');
end
ylim([0.3 0.40]);





